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We formulate a dynamical real space renormalization group approach to describe the time evo- 
lution of a random spin-| chain, or interacting fermions, initialized in a state with fixed particle 
positions. Within this approach we identify a many-body localized state of the chain as a dynamical 
infinite randomness fixed point. Near this fixed point our method becomes asymptotically exact, 
allowing analytic calculation of time dependent quantities. In particular we explain the striking uni- 
versal features in the growth of the entanglement seen in recent numerical simulations: unbounded 
logarithmic growth delayed by a time inversely proportional to the interaction strength. The particle 
number fluctuations by contrast exhibit much slower growth as log log t indicating blocked particle 
transport. Lack of true thermalization in the long time limit is attributed to an infinite set of ap- 
proximate integrals of motion revealed in the course of the RG flow, which become asymptotically 
exact conservation laws at the fixed point. Hence we identify the many-body localized state with 
an emergent generalized Gibbs ensemble. 



What is the effect of interactions on Anderson local- 
ization? One common wisdom is that any amount of in- 
teraction will give rise to collective excitations that could 
assist transport at non- vanishing temperature even if sin- 
gle particle states are all localized. But the belief, that 
there are no strict many-body insulators at T > 0, has 
been challenged by theoretical arguments, dating as far 
back as Anderson's original paper, which suggest a many- 
body localization transition marking a critical point in 
the transport properties of a closed quantum system [1, 2]. 
The idea has recently gained support from numerical 
studies [3-6]. Furthermore, simulations of one dimen- 
sional systems have revealed remarkably universal behav- 
ior of the dynamics in the putative many-body localized 
state[7-9]. For example the time evolution following a 
quench from a state with fixed particle positions shows 
blocked particle transport accompanied by unbounded 
logarithmic growth of the entanglement entropy. 

In this paper we formulate a real space renormalization 
group (RG) scheme that describes the time evolution of 
a random spin chain on multiple time scales. For certain 
initial conditions we can establish a many-body localized 
state as an infinite randomness fixed point of the dynam- 
ics near which the RG scheme is asymptotically exact. 
The results of this theory exhibit many of the universal 
features observed in the numerical simulations [7-9], al- 
though they are derived for a somewhat different model. 
As our starting point we consider the Hamiltonian of the 
random spin- 1/2 XXZ chain without local Zeeman fields: 

H = E I ( s ? s ^+i + s ^ s ^ + 2A ^+0 ■ (i) 

i 

The couplings Ji and anisotropy parameters Aj on sites 
i are random variables drawn from uncorrelated prob- 
ability distributions. The couplings may be positive or 
negative and we assume |Aj| < 1. The spins in (1) can be 
mapped using a Jordan- Wigner transformation to spin- 



less fermions with nearest neighbor interactions, subject 
to bond disorder. 

To study the propagation of information through the 
chain we investigate the time evolution of the system 
starting from a non entangled initial state, which for sim- 
plicity we take as an antiferromagnetic Neel state with 
spins pointing along the z-axis. We shall see that this 
choice of initial state greatly simplifies the scheme and 
allows us to obtain well controlled results for the dynam- 
ics at long times. 

RG scheme - Before proceeding with the details, let 
us outline the idea underlying the RG solution of the 
time evolution. As in the standard strong disorder RG 
scheme[10-12], we take advantage of the large local sep- 
aration of energy scales induced by the randomness to 
gradually eliminate degrees of freedom. However, instead 
of targeting the ground state of the system, our aim is 
the time evolution of the chain starting from a specified 
initial state. The dynamics at the shortest time scales 
are oscillations of frequency performed by the most 
strongly coupled pairs of spins on the chain, which are 
effectively decoupled on these time scales from their typ- 
ically much slower neighbors. On time scales longer than 
f2 _1 the rapid oscillations performed by the strong bonds 
can be eliminated using time dependent perturbation 
theory, which leads to renormalization of the coupling be- 
tween the slow spins. In this way we gain the essential in- 
formation about the dynamics of the chain at all scales. If 
the distribution of coupling constants flows to a wide dis- 
tribution at the dynamical fixed point, then the pertur- 
bative RG approach becomes increasingly well controlled, 
or even asymptotically exact if the system flows to infi- 
nite randomness [12]. Similar ideas have been applied to 
solve classical dynamics in certain disordered[13, 14] as 
well as clean[15] systems, but to our knowledge not to 
quantum dynamics. 

We now apply this scheme to the Model (1) with a 
staggered (Neel) initial state. The strong bond, with ex- 
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change coupling denoted by il, then always connects a 
pair of anti-aligned spins. At time scales of order SI -1 
the strong pair is effectively decoupled from the neigh- 
boring spins that are presumably connected to it by much 
weaker couplings. The only significant dynamics in this 
neighborhood of the chain is the rapid oscillation of the 
strong pair with frequency between the | tl ) and | ) 
states. 

The effective Hamiltonian for the dynamics at time 
scales much larger than is derived using time depen- 
dent perturbation theory in the coupling of the strong 
pair to the rest of the chain. Specifically, we move 
to the interaction picture with respect to the Hamil- 
tonian of the strong pair and compute the evolution 
of the density matrix p(t) — uj(t)\ip ) (ip \Ui(t) up 
to second order, while averaging over rapid oscillations 
of frequency fi. The resulting time evolution can be 
matched term by term to an effective time evolution 
exp(— iH c fft) \ tpo) (V'o I exp(£i?efft). Therefore this pro- 
cedure amounts to derivation of an effective Hamiltonian. 

At this order of perturbation theory the | tt ) and 
| J4- ) are not populated and therefore truncated from the 
Hilbert space. The retained states | ±> = 2~ 1 / 2 ( | tl) ± 
| I'f) ) of the strong pair can be taken as the f/i states 
of a new pseudo spin variable S n , which initially points 
along the positive or negative x-axis. The effective 
Hamiltonian derived in this procedure for our model is 
then given by: 
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Because [H e g, S*] = the time evolution can be com- 



puted separately for each eigenvalue ±\ of S£, using H^ s 
that docs not depend on the operator S n . The different 
evolution under H^ s together with the fact that the new 
spin starts in a superposition of f and J, leads to entangle- 
ment between the effective-spin on the strong bond and 
its two neighboring spins. Full entanglement is generated 
after a time £ en t = 2£!/( Jl JrAs), set by the difference 
in the exchange constant in H^ s . This process will be 
important later on for computing the evolution of the 
entanglement entropy. 

Apart from generating entanglement, the difference be- 
tween the evolution given f n or \, n is not crucial for the 
subsequent dynamics in the sense that they both lead 
to the same recipe for renormalization of coupling con- 
stants. H^g and have the same form as the original 
hamiltonian and we can directly read off the coupling 
generated between Sl and Sr, neighboring the strong 
bond to the left and right, upon decimation of that bond: 
J = JlJr/Q and |A| = |Al||Aa|/4, where we neglected 
the linear A correction to J. This approximation will be 
justified a posteriori by the fact that A flows to zero. The 
renormalization of the exchange coupling is then identical 
to that found in the random Heisenberg chain at T = 
and leads to the Random singlet phase[10-12]. Note also 
that we keep only the absolute value of the anisotropy. 
The sign will randomize in the course of the RG flow be- 
cause it depends on the state of S*. These RG steps are 
now iterated as we gradually approach lower frequency 
scales and longer times. 

The RG steps are iterated to produce a flow of the 
probability distributions with decreasing frequency cutoff 
fi starting from the microscopic cutoff fl . Using the 
scaling variables £ = In j and (3 = — In | A | , and T = 
ln(f2o/^) = ln(fioi) we obtain the following equation for 
the joint probability distribution P((,(3;T): 



dP 

or 



dP 



+ p(0; r) / dp L dp R d( L d( R 6(( -Cl~ (r)S(P -Pl-Pr- ln4)P(C L , H\ T)P(Cr, Pr\ H- 



(3) 



where p(C; r) = / dpP((, /?; T) is the distribution of (• 
Note that even if initially the variables C and j3 are inde- 
pendent a correlation builds up in the course of renormal- 
ization in the same way as it is generated in the ground 
state [12]. 

By integrating over f3 we obtain an equation for T), 

|£ = |?+p(0;r) J°° d( L dt R 6((-c L -CR)p(tL)p({R). 

(4) 

This equation is identical to the flow leading to the 
random-singlet ground state and it is solved by the same 



ansatz[12], p((;T) = a(T)e- a ( F K with 

a(r) = fTW (5) 

Of course the above solution includes only partial in- 
formation on the fixed point of the dynamics. Impor- 
tant information for calculation of physical quantities is 
held in the conditional average of the interaction vari- 
able (3 given a value of ( on the same bond, (3((, T) = 
J °° dPfiP((, p; T)/p(C; r). Proceeding as in Ref. [12], we 
derive the equation for this moment by multiplying Eq. 
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FIG. 1. Schematic illustration of remaining spins and clusters 
of decimated pairs in the renormalized chain at time t. 

(3) by f3 and then integrating over /3 
d r In (0(C) p) = d c In (0(0 p) + 2 ^£ <0( C ') (6) 

Note that we have neglected the In 4 in (3). This is jus- 
tified near the fixed point since the typical (3 flows to oo. 
After substituting the solution for p((, T) in Eq. (6) we 
find the solution 

^(C) = i(a r + l)^l + f ^ r ^ (7) 

where <f> = (1 + V / 5)/2 is the golden ratio and bo is deter- 
mined by the initial condition. 

Another important ingredient for calculation of phys- 
ical properties is the distance between remaining spins 
(i.e. length of decimated clusters) at time t. An il- 
lustration of such clusters is shown in Fig. 1. Since 
the RG flow is the same as in Rcf. [12] we similarly 
obtain L r = (a Q T + l) 2 = [a ln(fi *) + I] 2 Which be- 
haves as \n 2 (fl t) at long times. From this result we 
can immediately infer the decay of the staggered moment 
which is simply the fraction of remaining spins at time 
t: m s — 1/L(t) = 1/ (aoln(Qot) + l) 2 ■ It is interesting 
to contrast this behavior with the decay of the staggered 
moment in the analogous quench of a clean XX Z model. 
In that case Ref. [16] found oscillatory decay of the stag- 
gered moment (for A < 1) with an envelope that decays 
exponentially in time. 

We note that in the non-interacting case = the 
RG flow generalizes in a straight forward way to an ar- 
bitrary Ising initial states rather than a Neel state. If 
the strongest exchange coupling happens to be between 
spins with the same orientations, this pair will freeze in 
their state, while quantum fluctuations lead to an effec- 
tive coupling between the neighbors that is the same in 
absolute value as that obtained for anti-aligned spins. 
Hence the distribution of the absolute values of the ex- 
change couplings is governed by the same RG equations. 
In the interacting case, such aligned pairs lead to further 
complications that are left for later work. 
Results - We now use the solutions of the flow equa- 
tions to compute the time evolution of important phys- 
ical quantities. Specifically to gain information on par- 
ticle transport and on thermalization of the system we 
consider the growth of the total particle number fluctu- 
ation and of the entanglement entropy in a sub-system 
corresponding to half of the chain. 

To compute the particle number fluctuations we note 
that the total particle number (or S* oi ) within a deci- 



mated pair of sites is a conserved quantity in the RG 
scheme. Therefore the only contribution to the particle 
number fluctuation in one half of the system comes from 
decimated pairs that reside on different sides of the in- 
terface. Each oscillating pair that cuts the interface adds 
1/8 to the number fluctuation on time average. Comput- 
ing the total particle number fluctuation then amounts 
to counting the number of decimated bonds that cuts 
the interface. Proceeding exactly as in the random sin- 
glet phase[17] we have N p ps J r dT'a{T') = In (r + l/a ). 
Hence the particle number fluctuation grows extremely 
slowly as {SN 2 ) = (1/24) ln(ln(fi *)) at long times. In- 
terestingly this result is independent of the interaction 
strength A. 

We now turn to the growth of the entanglement en- 
tropy,between the two halves (A and B) of the system, 
S = — tr pa log 2 pa- Consider first the simpler (non- 
interacting) case Aj = 0, where no entanglement is gen- 
erated between a decimated pair and the rest of the sys- 
tem. Then, the only source of entanglement between 
the two halves of the system are decimated pairs that 
cut the interface, each contributes a time average of 
S p = 2-1/ In 2 ps 0.557. The growth of the entropy 
is therefore similar to that of the particle number fluctu- 
ation 

So(*)«Spi]n(In(not) + l/ao) (8) 

We can generalize this result (for A, = 0) to a quench 
from an arbitrary Ising state with a fraction q of anti- 
aligned neighbors. Because q is an invariant of the RG 
and aligned pairs do not contribute to the entropy the 
prefactor in (8) changes to qS p . 

We turn to our main focus, which is the interacting 
system with A / 0. In this case a pair decimated at 
time t\ will eventually get entangled with the neighbor- 
ing spins according to Eq. (2) after a characteristic time 
^ent(^i) = 20x / '{J\ Ai) • In particular from the time of 
the quench, entanglement will be generated by the in- 
teraction only after a delay time tdelay = 2S!o/(^oAo) = 
2(fio / Jo)0-/ Jq) j where J§ = JoAo is the typical value of 
the bare interaction energy. 

The interaction-generated entanglement entropy found 
at time t originates from entanglement of pairs elimi- 
nated at an earlier time t± = t — t crA , which corresponds 
to the scaling parameter T\ = lnfioii- To estimate this 
contribution to the entropy we recall that spins on the 
renormalized chain at time t\ are separated by clusters 
of length L(Ti) of spins decimated at even earlier times 
and are oscillating at higher frequencies. By the time t 
that a pair of spins decimated at t\ entangles with their 
neighbors, the spins inside the decimated clusters must 
also be entangled with each other. So, it is safe to assume 
that by the observation time t entanglement propagates 
to a distance L(Ti) giving rise to entanglement entropy 
S ps 0.5L(r x ) ps 0.5(a ri + l) 2 . The factor 0.5 stems 
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from the number of available degrees of freedom: the 
two states with aligned spins in each decimated pair re- 
main unpopulated and therefore do not contribute to the 
entropy. To write this as a function of the time t we use 
the relation between t and t\ 1 



t — ti + t cnt — t\ 1 



JlJrA! 



2VL\ 



JlJrAi 



(9) 



We now take the logarithm of both sides and replace 
the scaling variables by their appropriate average values 
C ->■ l/a(ri) and ft j3(( = 0;Ti). Here the correla- 
tions between the random variables came in to play. We 
needed the average of on the bonds with strongest J 
(C = 0), which is different than the global average of ft 

Using the solutions (5) and (7) for the typical values 
we find r = 3ri + ^(a ri + l)* + 2/a + In 2, This 
equation should be inverted to express Ti as a function of 
r and used to obtain S(T) = 0.5L(ri(T)). In general, the 
inversion cannot be done analytically, but we can easily 
find the limiting regimes. At sufficiently long times the 
term rf dominates the right hand side and then we have 
acTi = [&o(T - 2/ao - ln2)] 1 /* - 1. On the other hand 
at shorter times, when the linear term dominates, we get 
I\ = i(T- ^ - f -In 2). 

The crossover time t* separating the two regimes de- 
pends on the initial conditions through the coefficients of 
the terms Ti and Tf. If &o 3> ao, that is for stronger 
disorder in hopping than in the interactions, we have 
t* = tdeiay exp (6(3&o/ao)^/ a o)- In the opposite regime 
bo ao the term Tf dominates from the outset and 
t* = ^doiay We can now write an expression for the 
growth of the entanglement entropy valid in the limiting 
regimes 



S(t) 



1 /ln(t/t d elay) 



2 V ln(fi /J : 

1 ( In ftAdelay) 

ln(l/A ) 



+ 1 9(t- t delay )6(t* - t) 



2/4, 

+ I ! 6{t-t*) 



(10) 



Interestingly, Eq. (10) gives the unbounded logarithmic 
growth of the entanglement entropy seen in the numer- 
ical simulations and even the delay of this interaction 
induced growth by a time that scales as the inverse in- 
teraction strength [9]. It is important to note however 
that the numerically simulated model differs from ours 
by having disordered Zeeman coupling, which may ulti- 
mately drive it to a different fixed point [18]. Neverthe- 
less we can attempt a comparison even of the more subtle 
form of the logarithmic growth by crudely accounting for 
this difference. Through second order perturbation the- 
ory the Zeeman strong disorder generates predominantly 
strong randomness in the hopping, but not in the interac- 
tions. Hence the natural regime to compare the numerics 
to our model is a <C b where the first line of (10) holds 
at any reasonable time-scales (i.e. below the exponen- 
tially long time t*). Even within this regime Eq. (10) 



predicts a crossover between growth of the entanglement 
entropy as hit at times t < (fio / Jo) 2 tdeia,y to growth as 
In 2 t at longer times. The slope of the logarithmic growth 
1/ ln(f2g/Jo) is completely independent of the interaction 
strength A, consistent with the universality identified in 
Ref. [9]. 

To complete the comparison, we note that the growth 
of particle number fluctuations which was computed 
above is much slower (~ mint) than the growth of the 
entanglement entropy. Again this is consistent with the 
numerical results, which implies that particle transport 
is essentially blocked. Indeed the particle number fluctu- 
ations are in general only a lower bound of the entangle- 
ment entropy [19]. 

Having found that the entanglement entropy increases 
without bound it is natural to ask if this leads to ther- 
malization. To address this issue let us consider the sat- 
uration of the entanglement entropy in a finite system, 
or in a finite sub-system of length L s . Eq. (10) im- 
plies that the entropy will approach its maximal value 



^delay GXp 



ln(A )L 



0/2 



Does 



Soo after a time t sat 

the saturation value correspond to a state in thermal 
equilibrium? 

Provided we start from a symmetric distribution of A; 
such that (Aj) = 0, then the initial Neel state has zero 
mean energy, exactly in the middle of the many-body 
energy spectrum. If this state thermalized following the 
quench, the entanglement entropy would have to satu- 
rate to its infinite temperature value of L. But as we 
have pointed out above, the RG flow implies a satura- 
tion entropy that is at most half of the infinite tempera- 
ture value because half of the degrees of freedom remain 
frozen in the dynamics. This fact is embodied in the in- 
finite set of emergent integrals of motion I p — (SfS^p, 
which account for the fact that a pair of decimated spins, 
never flip their relative orientation within the perturba- 
tive RG scheme. Note that this remains true even for 
generalized initial states allowing for aligned neighboring 
spins. In the case of an initial Neel state we also have 
the particle number on oscillating pairs (or (Sf + 5 l |)p) as 
additional integrals of motion. These emergent conserva- 
tion laws become asymptotically exact for well separated 
pairs of spins decimated at long times as the perturbative 
RG scheme becomes asymptotically exact near the infi- 
nite randomness fixed point. We conclude that the long 
time steady state of the chain with non vanishing interac- 
tion is characterized by the generalized Gibbs ensemble 
(GGE)[20], which describes thermalization within a sub- 
space constrained by the values of the emergent integrals 
of motion I p . 

It is interesting to note that the long time steady 
state attained by the non interacting state Aj = is 
markedly different. The extremely slow increase of the 
entanglement entropy as lnlnt given by Eq. (8) to- 
gether with the relation between length and time scales 
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In fiot = r ~ v L imply saturation of the entanglement 
entropy to 5oo « In L. This result, as well as the lnlnf 
growth of the entropy matches with numerical results ob- 
tained for the random transverse field Ising chain[21] that 
can be similarly described by a model of noninteracting 
fermions. 

So far we have used the RG approach to establish and 
characterize a many-body localized state. It is interesting 
to examine the criterion for validity of the RG scheme, 
which may indicate a transition to a different, possibly 
thermalizing and delocalized, state. Such a criterion can 
be obtained from the first term in Eq. (2). The pertur- 
bation theory at the base of the renormalization step is 
valid while {J L J R /9) < (J)(l - A|). This can be ex- 
pressed in terms of the average values of parameters at 
the beginning of the flow. Using (J /Sl ) — a a /(ao + 1), 
we obtain the criterion eto < A^~ 2 — 1. Recall that increas- 
ing ao corresponds to decreasing disorder. It is tempting 
to speculate that for weaker disorder ao > A~ 2 — 1 the 
dynamics will give rise to a normal thermalizing state. 
However the above requirement may be too strict. Even 



if the criterion qq < A c 



1 is not satisfied initially, 



this may be corrected in later stages of the flow as A 
independently flows to smaller values. It is therefore also 
possible that the transition to a different state occurs 
only for (A 2 ) « 1 independent of the disorder strength. 

Conclusion - Using a real space RG scheme formulated 
in real time, we gave a dynamical description of a many- 
body localized state in a random spin chain, equivalent to 
interacting fermions with random hopping. Within this 
approach the localized state is characterized by a flow to 
an infinite randomness fixed point. Solution of the flow 
equations allows us to characterize this state in a rather 
detailed way. The results are found to be in excellent 
agreement with recent numerical simulations done on a 
similar, albeit not identical, model[9]. 

Particle localization is manifest in the extremely slow 
growth ~ lnlni of the particle number fluctuations in 
half the system that is seen both in the interacting and 
non-interacting systems. The entanglement entropy S re- 
veals a dramatic difference between the Anderson local- 
ized state of noninteracting fermions and the many-body 
localized state established when interactions are present. 
In the non interacting system S grows together with the 
particle number fluctuation as S(t) ~ log log t and satu- 
rates to a non-extensive value ~ InL in a finite system. 
Interactions lead to much faster growth of S(t) as log 2 ^ t 
at long times, but they take effect only after a delay 
time that scales as the inverse of the interaction strength 
^dciay ~ 1/J Z - Furthermore, the \og 2 ^ t behavior seen in 
the long time limit is preceded by logt growth up to an 
intermediate time scale ty m ~ idciay(^o/^o) 2 3> ^delay It 
is interesting to note that the growth of entanglement as 
In 2 ' ^ t exceeds the upper bound ~ In t proved for non- 
interacting Anderson localized chains [22]. 

The RG flow toward the infinite randomness fixed 



point has direct consequences on the equilibration in this 
system. In a sub-system of length L the entanglement 
entropy saturates to an extensive value Soo ~ L, which 
is however smaller than it would reach had the system at- 
tained true thermal equilibrium. We attribute the lack of 
thermalization to an infinite set of emergent integrals of 
motion, which become asymptotically exact conservation 
laws near the infinite randomness dynamical fixed point. 
The dynamics of the system can therefore be viewed as 
thermalization within a GGE characterized by the emer- 
gent set of conserved quantities, a possibility suggested 
in Ref. [23,24] . Here we demonstrated that such a GGE 
emerges in a non integrable random system as a dynami- 
cal fixed point of the renormalization group and captures 
the essence of a many-body localized state. The nature 
of the critical point marking the transition to the normal 
thermalizing state remains an interesting question for fu- 
ture study as are generalizations of our scheme to more 
generic disorder models and initial states. 
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